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INVERSE TRANSONIC AIRFOIL DESIGN METHODS 


INCLUDING BOUNDARY LAYER AND VISCOUS INTERACTION EFFECTS 

I. Introduction 

This report covers the period 1 August 1931 to 31 January 1982. The 
primary task during this reporting period was the extension of the grid 
embedment technique for TRANCES to transonic cases with viscous interaction. 
In addition, the study of the massive separation model (SKANFP) to medium 
speed cases was continued. 

II. Personnel 

The staff assigned to the project during this reporting period were: 
Leland A. Carlson, Principal Investigator 
August — Approximately h time 
Sept. -Nov. — Approximately 1/8 time 
Dec. -- Approximately 1/4 time 
Christopher Reed, Graduate Research Assistant 
August - December — H time 

The work assignments during this reporting period were: 

(a) Grid Embedment (Carlson and Reed) 

(b) Massive Separated Flow (Carlson) 

III. Embedded Grid Studies 

In the last progress report^, the second phase of the embedded grid 
studies was discussed. In that phase, the embedded grid scheme was redone to 
include upwind differencing and, as necessary, time-like damping. 
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It utilized a body fitted orthogonal curvilinear coordinate system, as shown 
on Figure 1, embedded in the overall cartesian system; and the appropriate 
full inviscid potential flow equation was solved using SLQR sweeping from 
F to B and F to C respectively. The size and discretization of this embedded 
grid is variable and determined by user selected input variables. This 
technique has been successfully applied to incompressible, slightly 
supercritical, and supercritical cases. However, in its second phase form, 
this approach did not include the effects of weak viscous interaction. 

During the present repr t '^iig period, this embedded scheme has been 
extended to include the effects of laminar-turbulent viscous interaction. 

As before, this new code uses for the etrtaedded grid an orthogonal body-fitted 
coordinate system in which one grid line coincides with the surface of the 
airfoil. At this grid line, the no normal flow boundary condition is enforced 
in the inviscid case. In the viscous case, an appropriate tangency boundary 
condition is imposed at the location of the displacement surface. 

The boundary condition imposed on the outer boundary of the en^edded grid 
requires that the embedded grid solution interface smoothly with the main grid 
solution. An initial embedded grid solution based upon the main grid solution 
is also needed to start the embedded iterative scheme. In order to satisfy 
these two conditions, a linear interpolation is used to obtain a perturbation 
potential f value at each embedded grid point from the main grid ^ values. 

The i values around the outer embedded grid are then held constant during the 
embedded relaxation process to satisfy the outer boundary condition. 

The incorporation of weak viscous effects into the analysis of the airfoil 
flowfield is accomplished in the current program by assuming that the inviscid 
streamlines follow a displacement surface having ordinates and slopes different 
from the actual airfoil. 
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Figure 1 - Ccnputational Grid System used in TRAND£S 





The flowfield is then solved Inviscidly around this displacement surface 
to obtain the pressure distribution around the airfoil. 


The ordinates of the displacement surface are obtained by adding a 
calculated boundary layer displacement thickness to the airfoil ordinates. 

In order to calculate the boundary layer displacement thickness the 
pressure distribution around the airfoil must be known. Therefore, an 
initial inviscid flowfield solution must be obtained before viscous interaction 
can be introduced. .Once viscous interaction has been introduced, the 
flowfield must be solved iteratively along with the boundary layer in order 
to obtain the steady state solution. The values of displacement thickness 
actually used to obtain the new displacement surface ordinates are updated 
from one iteration to the next using under-relaxation, i.e., 

C-- 


where w is a relaxation parameter and S’* is the displacement thickness 
calculated using the latest pressure distribution. 

To obtain the boundary layer solution from the airfoil pressure 
distribution, the flow is assumed to be initially laminar at the leading 
edge stagnation point. In this region, a compressible Thwaites method, 
which is an efficient one parameter intregral method^ is utilized to obtain 
the boundary layer properties; and, once it has been determined that the 
flow has transitioned from laminar to turbulent, the viscous layer is solved 
using the Nash-Macdonald method with smoothing. The location of the transition 
point is determined from a Granville type correlation based upon the difference 
between the local momentum thickness Reynold's number and the value at the 
laminar instability point combined with the pressure gradient history. 
Sometimes, at high angles of attack, laminar separation is predicted upstream 
of this transition point on the upper surface. 
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If this situation occurs, the local momentum thickness Reynold's number is 
compared to an empirical correlation in order to determine if the resultant 
laminar bubble is long or short. Normally, the bubble is of the short type. 

If, however, it is long, the present model is not applicable and the subsequent 
calculations will probably be in error. In either case, transition to turbulent 
flow is assumed at the next grid point, and the calculation is continued. 

As stated, the Nash-Macdonald method together with certain smoothing 
operations is used to compute boundary layer properties in the turbulent 
flow region. In order to determine the boundary layer displacement thickness, 
the momentum integral equation must be satisfied for the momentum thickness 

This equation is . 

if *{Hi l-m‘) S : r 

and it is solved in conjunction with the formulas of Nash and Macdonald for 
skin friction, , and the shape factor, H = . This equation is numerically 

integrated at each grid point along the airfoil surface, and the resulting 
set of displacement thickness are then smoothed. At the trailing edge, an 
extrapolation procedure is used to obtain the trailing edge displacement 
thickness. The smoothing used in this procedure reduces the rapid variations 
which sometimes occur in regions with large pressure gradients, such as near 
shock waves. It should be noted that the Nash-Macdonald method with smoothing 
and extrapolation yields an overall and trailing edge behavior that agrees 
with the observed effect of the boundary layer on pressure distribution and lift. 

The procedure for inclusion of viscous interaction in the embedded grid 
case is similar to that used in the main grid where input to the boundary layer 
calculation routine consists of the pressure distribution and airfoil coordinates. 
Since the embedded grid consists only of points near the leading edge of the 
airfoil, the embedded grid solution must be interfaced to the main grid solution 
for input to the boundary layer calculation routine. 
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Each time the displacement surface Is updated, the pressure distribution 
for the embedded grid region Is calculated and used In that region; and, the 
pressure distribution from the last main grid solution Is used for the remainder 
of the airfoil. This new pressure distribution Is then used to calculate the 
airfoil boundary layer characteristics; and, the boundary layer displacement 
surface Is updated. 

The original reason for considering a coupled solution of the main and 
embedded grids was to correctly model the effects of very weak shock Induced 
viscous effects. These effects occur with the airfoil at a fairly high angle 
of attack In a medium (0.3-0. 5) Mach number flow. At these conditions, a 
small shock w^/e may be located very near the airfoil leading edge. Although 
this shock wave would be too small to detect using the main grid, the embedded 
grid should be able to capture It. This resolution should be very significant 
since even a small shock wave can have major effects on the subsequent 
development of the airfoil boundary layer. A coupled main-embedded grid 
solution Is necessary, therefore, to determine the effects of such a shock 
wave over the entire airfoil surface. 

In the cases previously reported the solution of the embedded grid 
always occurred after the main grid solution had been completed. However, 
to make a coupled solution possible, some method had to be developed In 
which the embedded grid Influenced the main grid. Several possible approaches 
were Investigated, but these did not yield any significant Improvements when 
compared to the main grid only solution. The present coupled method, however, 
allows small changes In the viscous boundary layer to affect the main grid 
solution, which Is physically realistic. 

In the coupled solution scheme, the main grid solution Is started In the 
usual fashion. However, at the time the viscous boundary layer update Is 
usually performed, an embedded grid solution Is obtained Instead. 
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This procedure Is followed by the viscous boundary layer update routine 
modified to include the embedded grid solution. After the boundary layer 
characteristics of the composite airfoil solution have been determined, a 
spline curve fit is used to update the displacement surface ordinates. 

The perturbation potentials and the other inviscid flowfield parameters from 
the previous main grid solution are not altered. The main grid solution 
routine is then continued with the updated displacement surface. A flow 
chart of this procedure is shown in Figure 2. 

Through numerical studies conducted as part of this research, an 
iteration procedure has been determined which rapidly and accurately 
converged to the final flowfield solution. First, two coarse Cartesian 
grids are solved completely inviscidly. These grids are typically 13x7 
and 25x13, although they can be varied by the user. These two Cartesian 
grids provide an initial solution to the finer Cartesian grids. On the 
third (typically 49x25) and remaining main grids, fifty inviscid iterations 
are completed initially, at which point an embedded grid solution is obtained. 
However, it is limited to one hundred iterations. Then, the embedded grid 
solution and the latest main grid solution are used together to determine 
the boundary layer properties. As noted previously, a spline curve fit updates 
the displacement surface ordinates. Subsequently, the main grid solution 
scheme continues for twenty more iterations, followed by another embedded 
grid solution and a displacement surface update. 

To increase the rate of convergence during initial boundary layer 
calculations, the relaxation parameter, w, is nominally set initially to 0.50. 
When the change in at the trailing edge from one viscous update to the 
next drops below ten percent, w is reduced to 0.25. 
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Figure 2 - Coupled Main-Enrt)edded Grid Solution Procedure 
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Frequently, the 5 * at the trailing edge converges very rapidly. In 
order to take advantage of this fact and to reduce computing time, the 
frequency of viscous updating Is reduced to every fifty Iterations once 
the change In at the trailing edge Is less than 0.0001. 

In order to further promote a converged solution, a certain nund)er of 
Iterations are set aside at the end of the procedure In which no displacement 
surface update? are computed. For the finest main grid selected, no updates 
are calculated after eighty percent of the maximum number of Iterations has 
occurred for that grid. For the previous main grid no updates are made after 
350 Iterations. 

Once the main grid solution has been reached, a final Inviscid eir^dded 
grid solution Is conducted so that the calculated lift and drag coefficients 
have the greatest amount of accuracy available. 

Several test runs have been conducted to verify the coupled main-embedded 
grid viscous solution method. The test cases chosen were for a completely 
subcritical flowfleld, a slightly supercritical flowfleld, and a highly 
supercritical flowfleld. In each case. In order to minimize computer costs^ 
only a 49x25 main grid with 25 points on each surface of the airfoil was used. 

The completely subcritical test case consisted of a NACA 0012 airfoil at 
a zero degree angle of attack and a Mach number of 0.10. The embedded grid 
used 79 points on the airfoil and extended rearward to x/c ■ -0.213. The 
pressure distribution for this case Is shown in Figure 3 along with the pressure 
distribution for a similar case using the main grid only for solution. Notice 
that the Cp's were slightly lower for the coupled solution than for the main 
grid alone solution. The drag coefficient due to friction (C^^) of the two 
solutions were also In good agreement. The of the coupled solution was 
0.0043 while that of the main grid only solution was 0.0049. 
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The trailing edge displacement thickness, i * of the main grid only solution 
converged tr 0.0036, while the coupled solution case trailing edge^* converged 
to 0.0029. 

The slightly supercritical test case consisted of a NACA 0012 airfoil 
at an angle of attack of one degree and a Mach number of 0.72. The embedded 
grid used 85 points on the airfoil surface and extended rearward to x/c=0. 18070. 
The pressure distribution of this case is shown in Figure 4 again compared to 
a similar result using only a main grid solution. Notice that the Cp's obtained 
from the main grid only scheme were higher than those resulting from the 
coupled solution procedure. As in the preceding case, the of the coupled 
solution was of the sante order of magnitude but slightly lower than the 
of the main grid only solution. In this case, the of the coupled solution 
was 0.0046 while the of the main grid only solution was 0.0048. The S* 

at the trailing edge of the two solutions also converged to similar results, 
with the coupled solution i* being slightly higher at 0.0064. TheS* at the 
trailing edge for the main grid only solution was 0.0060. In both cases, a 
short separation bubble was formed at the upper surface transition location. 

The location at which this ocurred was in very good agreement between the two 
solutions. The main grid only solution predicted the bubble and transition 
at x/c = -0.28287, while the coupled solution method predicted this bubble 
would occur at -0.28466. 

The highly supercritical test case consisted of a NACA 0012 airfoil at 
a Mach number of 0.80 and an angle of attack of zero degrees. The embedded 
grid used 79 points on the airfoil and extended rearward to x/c = -0.213. 

The pressure distribution of this case is shown in Figure 5 along with a similar 
result using only a main grid solution. 
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Notice that the Cp's obtained from the main grid only sch«ne were almost 
identical to those resulting from the coupled solution procedure. Again, 
the of the coupled solution was of the same order of magnitude but 
slightly lower than the of the main grid only solution. In this case, 
the of the coupled solution was 0.0037 while the of the main grid 
only solution was also 0.0037. However, the slight changes in the Cp 
distribution near the leading edge greatly improved the total drag. In the 
main grid only case, the total was only 0.0091 but with the embedded grid 
it was more reasonable 0.0128. The i* at the trailing edge for the two 
solutions again converged to similar values, with the coupled solution being 
0.0045 while the main grid only solution converged to a ^ * of 0.0046. The 
location of the transition from laminar to turbulent point was again in 
excellent agreement between the two solutions. The main grid only solution 
and the coupled solution both predicted transition at x/c = 0.09803. This 
agreement verifies that the introduction of the embedded grid viscous approach 
and its accompanying fine resolution enhances without degrading the accuracy 
of the boundary layer solution. 

Based upon these results, the following conclusions can be stated: 

1) A body-fitted grid embedment technique applicable to inviscid transonic 
airfoil flowfield analysis has been developed and verified through 

a series of test cases. 

2) Test cases used to verify the inviscid transonic airfoil flowfield 
grid embedment analysis technique show that the accuracy of the 
solution has been increased by grid embedding. This enhancement of 
the solution has been especially true when small supercritical zones 
occur which cannot be adequately described using the main grid only. 
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3) Viscous interaction has been applied to the body-fitted grid 
embedment analysis technique using the Thwaite's and Nash-Macdonald 
techniques. To obtain accurate results, however, it has been determined 
that the main and embedded viscous interaction solutions must be 
obtained by some type of coupled technique. 

4) A technique for solving the viscous, transonic airfoil flowfield 
using a coupled main and embedded grid has been developed and 
verified through a series of test cases. 

5) The test cases used to verify the coupled m^in-embedded grid solution 
technique show that the leading edge resolution is significantly 
increased without degrading the boundary layer solution. 

6) Results show that the embedded grid region can contain the entire 
supercritical zone or only a portion of it and still produce accurate 
results. It is noted, however, that placement of the main-embedded 
grid interface at the shock wave tends to degrade the solution. 

During the next reporting period, an attempt will be made to put existing 
embedded grid code into a more user orientated form. A rough draft of a user's 
manual has already been written, and it is planned to finalize this manual 
and issue it either as a NACA document or as a Texas Engineering Experiment 
Station report. 

IV. Massive Separation Studies 

In the last progress report^, several results obtained with the SKANFP 

full potential program were presented. This program is designed to handle 

massive separated flow and high lift and uses the simplified Kuhn-Nielsen 

2 

turbulent boundary layer method as modified by Barnwell . 
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As previously noted, the resultant Cp values exhibited an undesirable and 
unrealistic "bump” In the vicinity of the separation point due to a mismatch 
between the unseparated and separated pressure distributions. 

During this reporting period, several techniques have been studied In 
an attempt to eliminate this feature. After extensive investigation, this 
problem has been traced to the flowfleld solution at the point Just before 
separation, primarily in the value of the displacement surface slope just 
before the separation point. This value was previously determined by 
smoothing the displacement thickness, adding the smoothed values to the 
original ordinates, and spline fitting the resultant surface. The spline 
fit was then used to determine the derivatives of the displacement surface. 
Unfortunately, thei*'s used in this procedure included values in the 
separation zone predicted by the boundary layer calculations which may not 
agree with the shape determined by the inverse procedure in the separated zone. 

It was thought that perhaps the smoothing process was contributing to 
the problem since it would permit points in the separation zone to influence 
upstream points. Thus, a run with no smoothing on the upper surface 
displacement thickness was tried. Unfortunately, the results were oscillatory 
and unstable, indicating that smoothing was important to the overall process. 

Subsequent investigation indicated that an important quality affecting 
the results at the point before separation is the upper surface displacement 
slope. This value is "normally" determined by the spline fit routine. 

However, an alternative approach would be to compute it from the values of the 
displacement surface ordinate, and the U and V velocity components at thit 
location. This approach has been tried, and resultant slopes are compared 
to the original values on Figure 6. As can be seen, the original results 
(SKAN64) exhibit significant change just before separation; while the new data 
(SKAN66) is much smoother. 
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Figure 6 — Comparison of Upper Surface Slopes 
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This new approach looks promising; and It, among others, will be Investigated 
In more detail In the future. 

While the present results appear reasonable, significant problems still 
exist when they are compared to experimental data. Figure 7 compares the 
SKAN66 pressure distribution with test data obtained In the Langley LTPT 
wind tunnel. It appears that the present SKANFP method Is predicting separation 
too far downstream. This possibllty will be Investigated, and methods of 
correcting It will be developed. 

V. Publications 

During the present reporting period, the following publication was Issued: 
Reed, Christopher L., "Grid Embedment as Applied to Viscous Transonic 
Airfoil Flowfield Analysis", M. Sc. Thesis, Texas A&M University, December 
1981. 
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